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NOMENCLATURE 


dQr  d\ 

defined  by  eq  177  and  187 

/p  ,£(o*q  ) 

defined  by  eq  179  and  72 

^2  r 

defined  by  eq  188  and  199 

L 

latent  heat  of  fusion  of  water,  334  J/g 

dg, 

defined  by  eq  200 

L(ocle,a0) 

defined  by  eq  73 

A 

a  small  negative  number 

L(ais,a0) 

defined  by  eq  117 

Aq,  A-[ 

positive  numbers  used  in  eq  211 

Lc/L+ 

4_ 

defined  by  eq  126  and  162 

h 

positive  number  defined  by  eq  212-214 

L  L 
j-'c  /m: 

defined  by  eq  135  and  185 

where  i  =  1, 2, 3 

Eq,  Ep 

defined  by  eq  52  and  172 

Bi 

Ith  constituent  of  the  mixture.  Subscripts 
i  =  1, 2,  and  3  are  used  to  denote  unfro¬ 

mi 

positive  number  defined  by  eq  216  and 
217 

zen  water,  ice  and  soil  minerals,  respec¬ 

Aij  ,A4| 

names  of  models 

tively 

n 

boundary  in  R0 

c 

heat  capacity  of  the  mixture  defined  by 

h 

velocity  of  n  =  dn/dt 

eq  0 

boundary  with  i = 0, 1  where  n0  denotes 

co 

defined  by  eq  30 

the  boundary  where  T  -  0  (°C)  and  rq 

Ci 

heat  capacity  of  the  fa  constituent 

the  interface  between  R2  and  a  frozen 

Cs 

defined  by  eq  69 

fringe 

d 

unit  of  time,  day 

P 

pressure  of  water 

di 

density  of  the  fa  constituent 

Pa 

applied  confining  pressure 

e0r  el 

defined  by  eq  88  and  94 

Pn 

P(n) 

e2<  ^3 

defined  by  eq  119  and  112 

Po 

P  (no) 

e4,e4,es 

defined  by  eq  114, 145,  and  197 

q 

heat  flux  in  the  mixture  by  conduction 

e6>  el 

defined  by  eq  198  and  208 

Q 

defined  by  eq  5 

£\s>  ^os 

defined  by  eq  203  and  209 

defined  by  eq  91 

El,  E2 

defined  by  eq  121  and  122 

r 

rate  of  frost  heave 

E3,  E4 

defined  by  eq  139  and  140 

R0 

unfrozen  part  of  the  soil 

/ 

mass  flux  of  water  in  R, 

Ri 

frozen  fringe 

fo 

fi 

mass  flux  of  water  in  R0 

Rio 

part  where  0  >  T  >  Ta 

mass  flux  of  the  fa  constituent  relative 

Ru 

part  where  Ta  >  T  >  T1 

to  that  of  soil  minerals  where  i  =  1, 2 

r2 

frozen  part  of  the  soil 

fs 

defined  by  eq  68 

Si,  Si 

defined  by  eq  55  and  70 

F 

defined  by  eq  220 

Smr  Sp 

defined  by  eq  53  and  85 

got  gl 

defined  by  eq  127  and  99 

Sp,Sp, 

defined  by  eq  133  and  134 

git  83 

defined  by  eq  100  and  144 

Spp ,  Sx 

defined  by  eq  174  and  173 

hi,  h2 

defined  by  eq  63  and  64 

Si 

defined  by  eq  196 

hj,  /14 

defined  by  eq  77  and  78 

S2 

defined  by  eq  22 

H,  Hyj 

defined  by  eq  158  and  159 

S3 

defined  by  eq  88 

k 

thermal  conductivity  of  the  mixture 

t 

time 

K 

thermal  conductivity  in  R0 

T 

temperature  of  the  mixture 

h 

thermal  conductivity  in  R2 

Ti 

T^) 

K0 

hydraulic  conductivity  in  the  unfrozen 
part  of  the  soil 

Ta ,  Tb 

temperature  at  the  top  and  the  bottom 
of  a  sample 

Ki 

empirical  function  defined  by  eq  36 

Tc ,  Tm 

defined  by  eq  131  and  40 

where  i  =  1, 2 

Tp ,  Ts 

defined  by  eq  166  and  68 

K20 

defined  by  eq  214 

Ta,Tx 

defined  by  eq  37  and  97 

v 


Mi 

velocity  of  the  2th  constituent  where 
2  =  1, 2,3 

V 

defined  by  eq  16 

V0 

constant  speed  of  n0 

w0 

defined  by  eq  44 

Wi  ,  w2 

defined  by  eq  194  and  195 

W 

defined  by  eq  115 

w0,wx 

defined  by  eq  152  and  192 

X 

spatial  coordinate 

y>Y 

defined  by  eq  89 

z 

defined  by  eq  17 

a{t) 

trajectory  in  Figure  3 

a0 

absolute  value  of  the  temperature 
gradient  at  n0 

aj 

absolute  value  of  the  limiting  tem¬ 
perature  gradient  as  £  approaches  n\ 
while  ^  is  in  R2 

Ctie/ 

defined  by  eq  72  and  73 

o c0C,  aip  defined  by  eq  126  and  166 

y  constant,  1.12  (MPa/°C) 

8  thickness  of  a  frozen  fringe 

50 ,  rj  defined  by  eq  25  and  34 

0i  volumetric  content  of  the  ith  constituent 

rate  of  supply  of  mass  of  the  2th  constitu¬ 
ent  per  unit  volume  of  the  mixture 

ji  ,  ji0  defined  by  eq  104  and  116 
A  function  defined  by  eq  31 

v  ,  Vi  defined  by  eq  19  and  89 
Vi  defined  by  eq  157 

%  coordinate  defined  by  eq  10 

7t0 ,  defined  by  eq  44 

pi  bulk  density  of  the  2th  constituent 

Pio  PiinPo 

a ,  a0  defined  by  eq  38  and  24 
<rx  ,  <rc  defined  by  eq  129  and  130 
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Existence  of  Traveling  Wave  Solutions  to  the  Problem  of 
Soil  Freezing  Described  by  a  Model  Called  Mi 

YOSHISUKE  NAKANO 


INTRODUCTION 

The  scientific  study  of  soil  freezing  and  ice  segregation  began  in  the  early  1900s.  By  the 
1930s  researchers  (Taber  1930,  Beskow  1935)  had  already  found  that  ice  segregation  and 
the  resultant  frost  heave  are  caused  not  only  by  freezing  of  in-situ  water,  but  also  by  freez¬ 
ing  of  water  transported  toward  a  freezing  front  from  the  unfrozen  part  of  the  soil.  The 
understanding  gained  in  the  1930s  was  largely  qualitative.  However,  the  transport  of 
water  was  already  identified  as  one  of  major  issues  in  the  study  of  soil  freezing.  The  prob¬ 
lem  has  attracted  the  attention  of  many  researchers  (see  Nakano  1991). 

The  main  constituents  of  saturated,  frozen,  and  fine-grained  soils  are  a  solid  porous 
matrix  of  soil  particles  and  ice,  and  water  in  the  liquid  phase  called  unfrozen  water.  The 
physical  properties  of  all  constituents  except  unfrozen  water  are  well  understood.  It  is  gen¬ 
erally  understood  that  the  transport  of  water  in  frozen  soils  is  mainly  caused  by  the  move¬ 
ment  of  unfrozen  water  and  that  unfrozen  water  exists  in  small  spaces  surrounded  with 
surfaces  of  soil  particles  and  ice.  Heaving  during  freezing  is  not  limited  to  water  in  soil 
systems.  It  occurs  with  benzene  or  nitrobenzene  in  soils  (Taber  1930),  water  in  various 
powder  materials  including  hydrophobic  carborundum  (Horiguchi  1977),  liquid  helium  in 
porous  glasses  (Hiroi  et  al.  1989),  water  in  hydrophobic  silicon-coated  glass  beads  (Sage 
and  Porebska  1993),  and  water  in  porous  rocks  (Miyata  et  al.  1994). 

The  dynamic  and  thermodynamic  properties  of  liquids  have  been  known  to  be  modi¬ 
fied  by  confinement  in  very  small  spaces,  such  as  porous  media,  cell  membranes,  etc.  The 
problem  of  confined  liquids  has  attracted  the  attention  of  researchers  in  many  disciplines 
in  recent  years  (Granick  1991).  The  maximum  size  of  confining  space  that  significantly 
modifies  the  property  of  liquid  evidently  depends  on  a  kind  of  liquid  and  its  confining 
solid.  For  instance,  in  the  case  of  thin  quartz  capillaries  with  sizes  of  the  order  of  a  micron, 
the  melting  point  of  ice  is  practically  the  same  as  that  of  bulk  ice  (Churaev  et  al.  1993). 
However,  in  much  smaller  capillaries  of  the  order  of  50  nm,  the  melting  point  changes. 

A  significant  modification  may  occur  in  the  dynamic  behavior  of  water  in  fine  porous 
media.  It  is  known  (Angell  1983)  that  the  temperature  dependence  of  the  self-diffusivity  of 
supercooled  water  can  be  described  by  a  critical  type  of  equation  with  a  singular  tempera¬ 
ture  just  below  the  homogeneous  nucleation  temperature.  Recently  Teixeira  et  al.  (1997) 
have  found  that  the  self-diffusivity  of  supercooled  water  confined  in  fine  porous  silica  cor¬ 
responds  to  that  of  supercooled  water  at  about  30°C  lower  temperature.  Pagliuca  et  al. 
(1987)  have  shown  empirically  that  the  gradients  of  pressure  and  temperature  are  two  in¬ 
dependent  driving  forces  of  water  flowing  through  various  noncharged,  fine  porous,  and 
either  hydrophilic  or  hydrophobic  membranes  with  pore  size  of  the  order  of  10-500  nm  at 
temperatures  above  the  bulk  melting  point. 


The  specific  surface  area  of  fine-grained  soils  is  on  the  order  of  20-200  m2  / g,  and  unfro¬ 
zen  water  is  known  to  exist  in  the  form  of  thin  films.  The  thickness  of  such  films  depends 
on  temperature  and  pressure,  and  is  estimated  on  the  order  of  10-100  nm  at  the  tempera¬ 
tures  around  -  0.1°C  under  atmospheric  conditions  (Ishizaki  et  al.  1994).  Unfrozen  water  in 
frozen  soils  is  one  special  case  of  a  wide  class  of  confined  liquids.  The  key  issue  underlying 
the  transport  of  unfrozen  water  is  deemed  to  be  the  dynamic  collective  behavior  of  water 
confined  to  small  spaces  in  frozen  soils,  which  depends  on  complex  solid-liquid  interac¬ 
tions. 

It  is  generally  accepted  that  a  thin  transitional  zone,  often  referred  to  as  the  frozen 
fringe,  exists  between  the  0°C  isotherm  (frost  front)  and  the  growing  surface  of  an  ice  layer. 
The  unfrozen  water  content  in  frozen  soils  under  equilibrium  conditions  is  routinely  meas¬ 
ured  by  nuclear  magnetic  resonance,  differential  scanning  calorimetry,  or  time  domain 
reflectometry.  However,  since  the  unfrozen  water  content  under  dynamic  conditions  is  dif¬ 
ficult  to  measure,  the  phase  composition  of  a  frozen  fringe  is  not  known.  Since  the  proper¬ 
ties  of  all  parts  except  the  frozen  fringe  are  understood,  the  dynamic  behavior  of  the  frozen 
fringe  has  been  one  of  the  major  subjects  in  the  study  of  soil  freezing  in  recent  years.  Since 
the  1960s,  many  mathematical  models  (Talamucci  1977,  Kay  and  Perfect  1988)  of  a  frozen 
fringe  have  been  proposed  on  the  basis  of  various  hypotheses.  With  the  widespread  use  of 
computers,  the  methods  of  numerical  analysis  became  very  popular.  However,  because  of 
the  paucity  of  basic  knowledge  and  the  complex  nature  of  the  problem,  these  numerical 
studies  have  not  been  effective  for  the  critical  evaluation  of  the  multiple  hypotheses  used. 

Around  1980,  two  important  semiempirical  models  of  soil  freezing  were  introduced  for 
engineering  applications:  the  segregation  potential  (SP)  model  (Konrad  and  Morgenstem 
1981)  and  the  Takashi  model  (Takashi  et  al.  1978).  Today  the  SP  model  is  widely  used  for 
engineering  in  Europe  and  North  America,  while  the  Takashi  model  is  the  standard  of 
engineering  design  in  Japan.  These  two  semiempirical  models  share  a  common  approach 
that  the  freezing  characteristics  of  a  given  soil  are  determined  empirically  under  certain 
quasi-steady  conditions,  where  a  frost  front  moves  with  a  constant  speed.  These  models 
also  share  a  common  weakness  of  requiring  one  or  more  empirically  determined  parame¬ 
ters.  These  are  known  to  depend  on  not  only  the  properties  of  a  given  soil  but  also  a  partic¬ 
ular  quasi-steady  condition  specified  by  given  thermal  and  hydraulic  fields.  The  empirical 
determination  of  such  dependence  is  elaborate  and  costly.  An  accurate  mathematical 
model  is  needed  that  provides  the  functional  dependence  of  parameters  on  pertinent  vari¬ 
ables  specifying  given  thermal  and  hydraulic  conditions  in  terms  of  well-defined  functions 
(or  parameters)  describing  the  properties  of  a  given  soil. 

As  the  1980s  were  ending,  there  were  many  mathematical  models  of  soil  freezing 
(Gilpin  1980,  O'Neill  and  Miller  1985,  Fowler  1989,  etc.),  but  they  all  suffer  from  the  com¬ 
mon  fault  of  little  or  no  experimental  verification.  Efforts  were  initiated  to  study  the  prob¬ 
lem  analytically  and  to  verify  the  hypotheses  used  in  the  analysis  by  comparing  the  prop¬ 
erty  and  the  behavior  of  solutions  with  empirical  findings.  Adopting  such  an  approach, 
Nakano  (1990)  introduced  a  mathematical  model  called  Mi.  This  model  was  shown 
(Nakano  and  Takeda  1991,  1994)  to  be  consistent  with  experimental  data  on  the  growth 
condition  of  an  ice  layer  without  overburden  load  (Takeda  and  Nakano  1990)  and  under 
load  (Takeda  and  Nakano  1993).  The  growth  process  of  final  ice  lenses  was  accurately 
described  by  Mi  (Nakano  1992,  Nakano  and  Takeda  1993).  Nakano  (1994b)  has  shown  that 
the  functional  dependence  of  SP  on  thermal  and  hydraulic  conditions  predicted  by  Mi  is 
consistent  with  empirical  findings  that  were  used  to  build  the  SP  model. 

According  to  the  Takashi  model  the  freezing  characteristics  of  a  given  soil  are  described 
by  two  empirical  formulas  that  specify  the  dependence  of  the  frost  heave  ratio  and  the 
water  intake  ratio  on  given  thermal  and  hydraulic  conditions.  Two  theoretical  equations 
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corresponding  to  Takashi's  formulas  are  derived  by  using  the  analytical  solution  of  quasi¬ 
steady  problems  (Nakano  1994a,  Nakano  and  Primicerio  1995).  Comparing  the  theoretical 
formulas  with  the  empirical  ones  for  Kanto  loam,  Nakano  (1996)  has  shown  that  Mi  is 
compatible  with  the  Takashi  model.  Studying  the  property  of  a  frozen  fringe  described  by 
the  Gilpin  model  (Gilpin  1980),  Nakano  (1997)  has  shown  that  the  Gilpin  model  is  essen¬ 
tially  one  special  case  of  Mi  and  that  it  is  too  restrictive  to  accurately  describe  the  behavior 
of  two  kinds  of  porous  media  studied.  Assuming  linear  temperature  profiles  in  both  frozen 
and  unfrozen  parts  and  neglecting  the  effect  of  changing  composition  in  the  frozen  fringe, 
Talamucci  has  solved  the  first  (Talamucci  1998a)  and  second  (Talamucci  1998b)  boundary 
value  problems  of  unsteady  soil  freezing  based  on  M^ 

In  this  work  the  problem  of  soil  freezing  is  studied  by  using  Mi.  We  will  show  that 
traveling  wave  solutions  to  the  problem  exist  and  describe  how  these  solutions  have  been 
used  for  the  empirical  verification  of  M^ 


BALANCE  EQUATIONS  OF  MASS  AND  HEAT 

We  will  consider  the  one-directional  freezing  of  soils.  Let  the  freezing  process  advance 
from  the  top  down  and  the  coordinate  x  be  positive  upward  with  its  origin  fixed  at  some 
point  in  the  unfrozen  part  of  the  soil.  We  will  treat  the  soil  as  a  mixture  of  water  in  the 
liquid  phase  B lf  ice  B2 ,  and  soil  minerals  B3.  The  bulk  density  of  Bx  is  denoted  by  p*  (x,f).  If 
dx  is  the  density  of  the  zth  constituent,  then  the  volumetric  content  0*  (x,t)  of  the  ith  constit¬ 
uent  is  given  as 


0i  =  pi/di. 

It  is  clear  that  the  sum  of  0j  should  be  unity,  namely: 


(1) 


01  +  02  +  03  =  1. 


(2) 


We  will  assume  that  the  density  of  each  constituent  remains  constant. 

We  will  assume  that  the  unfrozen  part  of  the  soil  is  kept  saturated  with  water  at  all  times 
by  using  an  appropriate  water  supply  device.  The  balance  of  mass  for  the  z>th  constituent  is 
given  as  (Nakano  1990) 


dtPi~  aJPiMi)  +  V  1,2,3 


dx 


(3) 


where  ux  {xrt)  is  the  velocity  of  the  z>th  constituent,  and  Aj(x,  t)  the  time  rate  of  supply  of  mass 
of  the  zth  constituent  per  unit  volume  of  the  mixture.  The  summation  convention  on  index 
i  is  not  in  force  here,  so  that  (p^)  represents  only  one  term.  Since  none  of  the  constituents  is 
involved  in  a  chemical  reaction,  we  have 


Ax  +  A2  =  0  and  A3  =  0.  (4) 

We  will  assume  that  the  constituents  are  locally  in  thermal  equilibrium  with  each  other 
and  that  the  heat  capacity  q  of  the  ith  constituent  and  the  latent  heat  of  fusion  of  water  L  do 
not  depend  on  the  temperature  T.  If  k  is  the  thermal  conductivity  of  the  mixture,  the  con¬ 
ductive  heat  flux  q  (x,t)  in  the  mixture  is  assumed  to  be  given  as 
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Using  eq  5,  we  will  obtain  the  balance  equation  of  heat  for  the  mixture  (Nakano  1990) 
given  as 

—q  =  L[X2  +  z)  (6) 

where  z  (x,t)  is  defined  as 


Lz  =  -c^+ (cj  -  c2)T  X2  -  IPiMiCi 3T 


C  =  CD  +C,P  +cp,. 
lrl  2r2  3r3 


dx 


(7) 

(8) 


We  will  consider  a  special  case  in  which  a  frost  front  x  =  n()  ( t )  moves  with  a  constant 
speed,  namely 


~n  (f): 
dt  °v  ' 


-n  =7  >0. 

O  O 


(9) 


Hereafter  we  will  exclude  the  case  of  negative  Vo  where  melting  occurs.  We  will  introduce 
a  new  independent  variable  ^  defined  as 


F=x- it  t-n  (0). 
^  o  o'  ' 


(10) 


For  the  sake  of  convenience  we  will  define  new  dependent  variables/i  (^)  and/2  (£,)  as 


a = pi(«i  -  »3)  (ii) 

f2  =  p 2(u2  -  u3).  (12) 

Therefore,  ( i  =  1,2)  is  the  mass  flux  of  either  B1  or  B2  relative  to  the  mass  flux  of  soil 
particles.  Using  eq  10, 11  and  12,  we  reduce  eq  3  to 

(13) 

(,4) 

(p3V)'  =  0  (15) 

where  primes  denote  differentiation  with  respect  to  £  and  V  (£)  is  defined  as 

V  =  u3-n0.  (16) 

Similarly  we  will  reduce  eq  6  and  7  to 

q'  =  -(kT')’  =l(X2+z)  (17) 

Lz  =  -(c^  +  c2/2  +  cV)T  +  (Cl  -  c2)X2T.  (18) 

QUASI-STEADY  PROBLEM 


A  freezing  soil  may  be  considered  to  consist  of  three  parts:  the  unfrozen  part  RQ,  the 
frozen  fringe  Rx  and  the  frozen  part  R2,  as  shown  in  Figure  1.  We  will  also  make  seven 
assumptions: 
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1.  The  dry  density  of  R0  remains  constant, 

2.  The  composition  is  continuous  at  n0, 

3.  The  pressure  P  of  water  at  n  remains  con¬ 
stant  at  Pn , 

4.  f2  vanishes  in  Ri  and  R2  unless  p3  vanish¬ 
es, 

5.  The  Aux/l  is  negligibly  small  in  R2, 

6.  Sensible  heat  terms  are  negligible  in  com¬ 
parison  with  latent  heat  terms, 

7.  px  is  given  in  Ri  and  R2  as 

pi=p3VC0-  (19) 

The  bulk  density  pi  under  equilibrium  condi¬ 
tions  is  known  to  be  given  by  eq  19  where  v  (T) 
is  an  empirically  determined  and  increasing 
function  of  T.  Hence,  the  assumption  7  implies 
that  p!  under  dynamic  conditions  is  also  given 
by  the  same  form  as  eq  19.  We  will  assume  that 
v  (T)  has  a  continuous  first  derivative. 

We  will  seek  a  traveling  wave  solution  to  the  problem  in  which  the  boundaries  n  (f),  n0 
(t)  and  rii  ( t )  move  with  the  same  constant  speed  V0,  namely 

V  =- n--n  =-n.  (20) 

From  a  physical  point  of  view,  maintaining  a  constant  pressure  Pn  is  difficult  at  the  moving 
boundary  n  (f).  However,  a  solution  obtained  under  such  an  idealized  condition  is  quite 
useful  for  applications  (Nakano  and  Primicerio  1995).  If  such  a  solution  exists,  it  must  sat¬ 
isfy  eq  13- 15,  and  eq  17. 

From  eq  13, 14,  and  15,  we  find  that  the  flux  of  water  fi  (£,)  is  given  in  i\  j  (Nakano  1994a) 
as 

A  =fo  +  s2  (p10  -  Pgo^o  - ' h  (V-  V0),  0  <  l <  5  (21) 

where  pio  and  p3o  are  the  constant  bulk  densities  of  Bi  and  B3  in  R0,  respectively,  and 
8  =  ri]  -  n0,  and/0  is  the  constant  flux  of  water  in  R0.  s2  is  defined  as 

s2=l-d~1d2.  (22) 

Neglecting  the  gravitational  effects  and  using  Darcy's  law,  the  flux  of  water/0  in  R0  is  given 
as 


T a 


o°c  • 


Figure  L  Quasi-steady  freezing  of  soil 


/0=Koco5o1  (23) 

where  K0  is  the  hydraulic  conductivity  of  R^.  a0  and  S0  are  defined  as 

°a=Pn-Po'  Po=P(°)  (24) 

80=M0-n.  (25) 

The  boundary  n\  is  a  free  boundary.  The  composition  may  be  discontinuous  at  and 
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the  limiting  value  p3  (8+)  of  p3  as  £  approaches  8,  while  £  is  in  R2  is  given  (Nakano  1994a)  as 

p3(8+)  =  P3oVoW8+)  +  Vor1  (26) 

where  r  (8+)  is  the  rate  of  heave  at  8+  given  as 

r(8+)  =  d~lf{ •  +  d21s2[p10  -  p3ov{T(8)}F0.  (27) 

The  heat  flux  is  discontinuous  at  ?i;  and  the  jump  condition  is  given  as 

9(8+)  =  <7(8-)  +/i(8-)[L  +  {ci  -  c2)Ti ]  (28) 

where  g(8  -)  and/j(8  -)  are  limiting  values  of  cj  and/i,  respectively,  as  £,  approaches  8,  while 
£  is  in  Ri  and  Tx  is  T  (8). 

We  will  reduce  eq  17  and  eq  18  to  a  simpler  form.  Using  eq  13, 14,  and  15,  we  obtain 


cV  =  cqVa-(c1-c2)A  +  c1(f0-f1)  (29) 

where 

co  =  clPlO  +  c3p30  (30) 

A®  =  J  \d%,  ^>0.  (31) 

0 

Using  eq  29,  neglecting  a  sensible  heat  term,  and  integrating  eq  17,  we  obtain 

—k  T'  =  k0a0  +  LA,  0<^<8  (32) 

where  k0  and  k  are  the  thermal  conductivities  of  R0  and  Ri,  respectively,  and  a0  >  0  is  the 
absolute  value  of  the  temperature  gradient  at  £  =  0.  Using  eq  32  and  neglecting  sensible 
heat  terms,  we  will  reduce  eq  28  to 

*!<*!  -  koao  =  [^(S-)  +  A(8-)]l  (33) 

where  k\  is  the  thermal  conductivity  of  R2  and  cxi  is  -T '  (8+). 


Using  the  principle  of  mass  and  heat  conservation,  we  have  derived  equations  that  must 
be  satisfied  by  a  traveling  wave  solution  of  soil  freezing.  Clearly  these  equations  are  not 
sufficient  to  solve  the  problem.  We  need  a  model  of  a  frozen  fringe  that  specifies  fi  (i;)  and  A 

MODEL  STUDY 

A  model  of  a  frozen  fringe  called  Mi  was  introduced  by  Nakano  (1990)  to  explain  empir¬ 
ical  findings  on  the  growth  condition  of  an  ice  layer  in  freezing  soils.  The  model  has  been 
modified  as  its  empirical  evaluation  has  progressed  (Takeda  and  Nakano  1990,  Nakano 
and  Takeda  1991,  Takeda  and  Nakano  1993,  Nakano  and  Takeda  1994).  The  latest  version 
assumes  the  validity  of  equations  in  Ri  given  as 

k  -  constant,  k  /  k°  =  r|  >  1  (34) 
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P3=P30'  pl=  ^30^10'  p1(°+)  =  p10 

(35) 

/■/1=-K1P'-K27\  1^(0+)  =  ^ 

(36) 

K2(T)/K1(T)  =  y  for  0>T>T=-a/y 

(37) 

P(S-)  =  P  =a+P  ,  <j>0 

(38) 

P'(S-)  >  0,  Vo>0  and  V0P'(g~)  =  0 

(39) 

K2(T)  =  0,  K2(T)/K1(T)  =  0  for  T<Tm<Tc 

(40) 

where  7  is  a  constant  (1.12  MPa/°C),  P  (£)  is  the  pressure  of  water,  Pa  is  the  applied  confin¬ 
ing  pressure  (uniaxial  stress),  a  is  the  effective  confining  pressure,  and  K{  (i  =  1,2)  is  the 
transport  property  of  a  given  soil  that  generally  depends  on  the  temperature  and  the  com¬ 
position  of  the  soil.  Since  p3  is  a  constant,  we  will  assume  that  Kj  is  an  increasing  function  of 
T  alone.  This  assumption  implies  the  homogeneity  of  soils  in  a  microscopic  scale  that  corre¬ 
sponds  to  the  thickness  of  the  frozen  fringe,  which  is  clearly  an  approximation.  We  will 
assume  that  (T)  has  a  continuous  first  derivative.  Because  of  eq  37,  K2  (T)  may  be  discon¬ 
tinuous  at  T  =  Ta.  We  will  assume  that  the  first  derivative  of  K2  is  continuous  except  at 
T  =  T0.  It  is  known  that  the  mobility  of  unfrozen  water  tends  to  diminish  as  T  decreases.  We 
will  assume  that  there  exists  a  negative  number  Tm  <  Ta  such  that  eq  40  holds  true  and  that 
K2(T)>0  and  K2  (T)/  Ki(T)>  0  for  T  >  Tm.  According  to  Ml7/is  given  by  eq  36  in  i?i  while 
Darcy's  law  holds  true  in  R0 .  Hence,/ and  P  are  continuous  but  P'  maybe  discontinuous  at 

H0. 

The  Mi  model  is  a  generalization  of  somewhat  simpler  but  more  restrictive  models.  Mi 
(Derjaguin  and  Churaev  1978,  Ratkje  et  al.  1982,  Horiguchi  1987),  in  which  the  ratio  K2/Ki 
is  equal  to  y  regardless  of  T.  In  Mi  the  coupling  mechanism  for  mass  and  heat  transport  is 
based  on  irreversible  thermodynamics  in  which  local  equilibrium  is  assumed  under  a  tem¬ 
perature  gradient  (Ratkje  and  Hafskjold  1996).  In  Mi  local  equilibrium  holds  in  the  part  Ri0 
where  Ta  <  T  <  0,  but  not  in  the  part  Rn  where  T  (8)  <  T  <  Ta  (Fig.  1).  This  generalization  is 
needed  because  Mi  is  too  restrictive  to  accurately  describe  the  behavior  of  porous  media 
(Nakano  1997).  Equation  39,  often  referred  to  as  the  Signorini-type  free  boundary  condi¬ 
tion  (Friedmann  and  Jiang  1984),  is  needed  for  the  uniqueness  proof  of  solutions  when  VQ 
is  positive.  It  is  not  certain  that  such  a  condition  holds  true  because  of  the  paucity  of  exper¬ 
imental  data.  In  addition  to  the  above  equations,  we  will  assume  that  the  thermal  conduc¬ 
tivities  k0  and  are  given  constants  for  the  sake  of  simplicity. 

When  eq  35  holds  true,  u3  vanishes  and  eq  21  is  reduced  to 

/©  =/o  +  s2(PlO  “  P30v)^o/  0<^<8.  (41) 

The  A  (£)  is  given  as 

A(^)  =  d1-Id2(p10-P3ov)vo.  (42) 

According  to  M\  the  properties  of  a  given  soil  are  described  by  three  empirically  deter¬ 
mined  functions  of  T :  K\,  K2  and  vthat  are  assumed  to  be  functions  of  T  alone  for  T  <  0°C. 
The  hydraulic  field  is  specified  by  Pn,  80  and  Pa  while  the  thermal  field  is  specified  by  a0 
and  Oj.  Our  problem  is  to  find  constants  yo  >  0,  5  >  0  and  functions/©,  T  ©  <  0,  P  ©  so 
that  the  following  equations  (PI  through  P 7)  are  satisfied: 
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From  eq  41  we  have 


m  -  4 + *2  [pI0  -  P„v{TtO}] V  ,  0  <  5  <  8.  (PI) 

From  eq  36  we  have 

m=-K1{m)}p\^)-K2{m)}nk) ,  o<^<s.  (P2) 

From  eq  32  and  49  we  have 

kT%)  =  -koao-d?d2[pw-p30v{T&)}]LVo,  0<$<5.  (P3) 

From  eq  33  and  42  we  have 

fcxot!  -  fc0a0  =  L/0  +  [p10  -  p30v{T(8)}]LVro.  (P4) 

Boimdary  conditions  are  given  as 

T  (0)  =  0  (P5) 

P  (8  -)  =  Pa  (P6) 

P'(8~)  >  0,  P0  >0  and  VoP'(8-)  =  0.  (P7) 

We  will  rewrite  eq  P3  as 

r\T'  =  -jij  +  7t0v  (43) 

where  n0  and  nj  are  given  as 

7I0=df1d2fco1P3oPPo.  Jll=a0+7t0W0,  w0  =pio  /  p30-  C44) 

Since  v  (T)  <  zv0  for  T  <  0°C  and  a0  >  0,  T(t,)  is  strictly  negative.  Hence  the  function  T  (£)  is 
invertible  for  5  >  \  >  0.  Integrating  eq  43  by  using  eq  P5,  we  obtain 

m)  =  -(n1/^-(Kjn1)jT0v[l-(no/n1)v  ]~'dT.  (45) 

Integrating  eq  P2,  we  obtain 

P[£(T)]  -  Pn  +  (80  /  fC0)/0  =  -JoT  (K2  /  Ki)dT  -JoT  /[^(T)](K1T')-1  dT.  (46) 

Setting  T  =  Tj  in  (46)  and  using  eq  P6,  we  obtain 

+  (80  /  Ko  )/o  =  - (*2  /  Ki )dT  + /^DftKiT')-1  dT.  (47) 


Equation  47  provides  the  functional  dependence  of  Tj  on  a0  and  ai  that  specifies  a  given 
thermal  condition  as  well  as  on  80  and  ct  that  specifies  a  given  hydraulic  condition  in  terms 
of  functions  and  parameters,  such  as  K2,  K2,  and  v,  etc.,  describing  the  properties  of  a  given 
soil. 
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GROWTH  OF  ICE  LAYERS 

We  will  seek  solutions  in  which  V0  -  0.  In  this  case  eq  PI  through  P4  are  reduced  to 


M)=fo  ,  0<S<5  (48) 

m  =  -Kt{T(5)}P'(£)  +  K2{m}(aQ  /  ti)  ,  0  <£,  <  8  (49) 

TW  =  -(< a0/H)  (50) 

fcicci —  k0oc0  =  Lf0.  (51) 


The  left-hand  side  of  eq  51  is  the  rate  of  heat  removal  from  the  frozen  fringe  that  must  be 
positive  during  soil  freezing.  Hence,/0  >  0.  We  will  consider  a  quadrant  S  =  [(on,  a0):  on  >  0, 
a0  >  0],  where  we  draw  a  straight  line  Le  starting  from  the  origin  (Fig.  2a)  defined  as 


Le  ={(a!,a0):a0  =  (k1/fc0)a1}. 


(52) 


a.  Stable  growth  region  of  an  ice  layer. 


cn — *• 

c.  Region  of  frost  penetration  when  a>ac 
Figure  2.  Temperature  gradients  a-j  and  aQ. 


9 


It  follows  from  eq  51  that/0  vanishes  on  Le.The  line  Le  divides  S  into  two  regions,  and  we 
will  denote  one  of  them  by  Sm  defined  as 


Sm  =  {(ai'ao):a0  >  (h  /  fc0)ai}.  (53) 

Therefore,  melting  takes  place  in  Sm.  We  will  exclude  Sm  from  our  discussion  hereafter.  For 
a  special  case  where  a0  =  0  and/0  >  0,  from  eq  49  we  find 


/o=-Ki(ri)P'(5-).  (54) 

It  follows  from  eq  54,  P 7  and  51  that  cq  and/0  also  vanish.  Now  we  will  seek  solutions  with 
V0  =  0  in  the  region  Sf  defined  as 

Sf  =S-(Sm+Le)  =  {(a1,ao):0<ao  <(fc1/fc0)a1}.  (55) 

Suppose  that  such  solutions  exist  in  Sf.  Then,  eq  26  and  27  are  reduced  to 

p3(5+)  =  0  (56) 

r(5+)  =  d~1fo>0.  (57) 

It  follows  from  eq  56  and  57  that  an  ice  layer  grows  in  the  solutions 
When  V0  =  0  in  Sf ,  eq  47  is  reduced  to 

°  +  (50 /K0]fo  =  £°  (K2  /KOdT-ft/ a0)/0£  (K.T'dT.  (58) 

We  will  write  eq  58  as 

(80  /  K0)  +  (Tl  /  a0)j£  (Ki)-idT]/0  =  (K2  /  Ka)  -  a.  (59) 

Suppose  that  >  Ta ,  using  eq  37,  we  obtain 

£(K2/K1)dT-o  =  y(Ta-Tl)<0,  if  Ta>Ta.  (60) 


It  follows  from  eq  59  and  eq  60  that/0  <  0.  Therefore,  if  there  exists  7\  such  that/0  >  0  and  eq 
58  holds  true,  then  T1  must  be  less  than  TCT.  Our  next  aim  is  to  find  points  (cq,  aQ)  in  Sf  with 
Tf  <  Ta  and  F0  =  0.  Hereafter,  we  will  assume  that  <j  >  0  and  80  >  0  are  given  constants. 

By  eq  P7  P'  (5  -)  is  nonnegative  if  V0  =  0.  We  will  begin  our  study  with  a  special  case  in 
which  P'  (5  -)  vanishes. 

Proposition  1 

If  P'  (5  -)  vanishes  for  a  given  a0,  then  there  exists  a  unique  Tj  such  that  Tm  <  Tf  <  Ta  and 
that  eq  58  holds  true. 

Proof 

When  P'  (8  -)  vanishes,  then/0  is  given  as 

/o  =  K2(T1)(a0/n).  (61) 
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Using  eq  61,  we  will  write  eq  58  as 


^i(^i)  -  ^2(^1)  (62) 

where  \  and  h2  are  defined  as 

fciCTi)  =  CT  +  (50/K0)K2(T1)(a0/n)  (63) 

h2 fli)  =  J°  [K2 (s)  -  K2(Ta)][l  /  Kjfe)]*.  (64) 

Since/0  >  0,  then  Tx  <  Ta.  Using  eq  37,  we  will  reduce  eq  64  to 

W  =  a  +  l^ K(s) /  K1(s)]ds-K2(ri)jJ  [l /  Kj (s)]ds.  (65) 

Differentiating  fej  (T{)  and  h2  (T-i)  with  respect  to  Ti,  we  obtain 

MT1)  =  (S0/K0)k2(T1)(a0/r\)  (66) 

h2(T1)  =  -K2(\)l0Til1/KMds  (67) 


where  a  dot  denotes  differentiation  with  respect  to  TV  Therefore,  from  eq  63  and  65,  hx  (Tm) 

<  h2(Tm)  and  /ii(Tff-)  >  h2(Ta-).  From  eq  66  and  67  we  find  that  h\(Ti)>  0  and  h2(Ti)  <  0  for 
Tm  <  Tj  <  Ta.  Also  hi(T{)  and  h2(T{)  are  continuous  for  Tj  <  T0.  Therefore,  there  exists  a 
unique  T1  such  that  Ti  <  Ta  and  that  eq  44  holds  true.  □ 

We  will  denote  the  unique  Ti  of  Prop.  1  for  a  given  a0  by  Ts  (a0).  It  is  easy  to  see  that  the 
function  Ts  (a0)  is  continuous  for  <  T0.  If  we  denote/0  by/s,  when  T1  =  TS,  then/s  is  given 
as 


fs  =  K2(Ts)(  Oo/n). 


(68) 


Substituting/0  in  eq  51  with/s,  we  will  define  a  curve  Cs  in  Sf  (Fig.  2a)  as 

Cs  =  {(ai,  a0):a0  =  kx{k0  +  -n-iJLiC2{Ts(cx0)}]-ict1 } .  (69) 

We  also  define  the  region  Si  bounded  by  Le  and  Cs  as 

Si  =|(ai,a0):(fci  /  k0)ax  >  a0  >  kx[k0  +Tr1LK2{Ts(a0))]_1a1}.  (70) 

From  eq  62,  Ts  depends  on  a0,  50,  and  a  for  a  given  soil.  We  will  show  the  nature  of  such 
dependence  below. 

Proposition  2 

The  solution  Ts  of  eq  62  is  a  decreasing  function  of  a0,  S0  r  a080  and  a. 

Proof 

Differentiating  eq  62  with  respect  to  a0,  we  obtain 

K2 (Ts)g-{5o(pK0 )-1a0  +  J®  [l /  Kx (s)]*}  =  -80(qK0 )-i K2 (Ts ).  (71) 
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It  follows  from  eq  71  that  Ts  is  a  decreasing  function  of  a0.  Similarly  it  is  easy  to  find  that  Ts 
is  a  decreasing  function  of  80 ,  a080  and  a.  □ 

Next  we  will  study  the  region  Sj.  For  a  given  a0,  we  will  consider  a  segment  ^(cx0)  of  a 
straight  line  L(aie,a0)  (Fig.  2a)  defined  as 

t(a0 )  =  {(oci,a0)  :  a0  =  a0  and  aie  <  ai  <  ajs}  (72) 

L(aie,a0)  =  {(a1,a0) :  a0  =a0,a!  >ale}  (73) 

where  (aie,  a0)  e  Le  and  (ais,  a0)  e  Cs.  The  flux/0  (aj)  on  i(a0)  is  given  by  eq  51  as 


f0(a1)  =  (k1al-k0a0)/L.  (74) 

The  flux  depends  linearly  on  air  /0( ai)  >  0  on  t(a0),  /0(ale)  =  0  and  /0( ais)  =/s  where/s 
is  given  as 

/s  =  K2{Ts(a0)}(a0  /r\).  (75) 

The  P'  (8  -)  vanishes  at  the  point  (ais/a0).  Suppose  that  P'  (8  -)  vanishes  at  some  point  on 
£(a0),  then  there  is  no  solution  of  eq  58  at  that  point  by  Prop.  1.  We  will  seek  solutions  on 
4«o)  under  the  condition  of  P'  (8  -)  >  0. 

Proposition  3 

For  a  given  a0,  there  exists  a  unique  on  £(a0)  such  that  TS<T1<  Ta  and  eq  58  holds 
true  if  P'  (8  -)  >  0. 

Proof 

For  a  given  point  (cy  ,cx0)  on  £(a0)  we  will  write  eq  58  as 


h3(a1)  =  h4(T1,a1)  (76) 

where  h3  and  /i4  are  defined  as 

^(ai)  =  ^+(8o/K0)/0(a1)  (77) 

h (Ti ,  (Xi)  =  J®  [K2  (s)  /  Kj (s)]ds  -  (n  /  Oo)/0 (cd)J®  [1  /  Ki(s)]*.  (78) 

Differentiating  /z4  (T4,  a4)  with  respect  to  T1(  we  obtain 

ft4(Ti/ai)  =  -(Ti/a0)P/(6-).  (79) 

Since  P'  (5  -)  >  0,  oq)  is  a  decreasing  function  of  Ty  When  approaches  Ts(cx0) ,  we 
will  evaluate  ^4.  Since  Ts(a0)  is  the  solution  of  eq  62,  from  eq  63  and  64  we  obtain 

a  +  (Sc  /K0)/s  =  J®  [K2(s)-X2(rB)l|l/ X!(s)lds.  (80) 

Using  eq  80,  from  eq  78  we  obtain 

fe4(Ts,a1)  =  o+(80/K0)/s+(i1/a0)K1(Ts)P'(8-)^0  [l/JC^s)]*.  (81) 


It  follows  from  eq  80  and  81  that  h^(T$/  oti)  >  /^(oti).  Also  when  T\  approaches  Ta,  we  have 
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(82) 


h (To  ,  a1)  =  c-(ri/ a0)/0  f  [1  / (s)]ds < h3 (ai)] . 

Therefore,  there  exists  a  unique  7)  such  that  Ts  <  <  Ta  and  eq  58  holds  true.  □ 

Differentiating  eq  76  with  respect  to  oq,  we  obtain 

(L  /  =  -a050 (K0-n)-i  -  J°  [1  /  X!(s)]dS .  (83) 

From  eq  83  we  find  on  ^(a0) 

dT\  n  /n,4\ 

^-^cO,  and  L->-cx>  as  ^  ^  als.  (o4) 

aCXi  dai 

As  oti  increases  from  aie  to  ais  on  ^(a0),  the  flux/0  increases  from  zero  to/s,  while  Ti 
decreases  from  Ta  to  Ts.  In  Proposition  3  o:0  is  an  arbitrary  positive  number.  Hence,  we  may 
conclude  that  an  ice  layer  grows  in  the  region  Sj  and  on  Cs.  Below  we  will  study  the  region 
Sp  defined  as 

Sp  =  |(ai,a0):fci[fc0  +ir1LK2(Ts(a0)}]-1ai  > a0  >  o}.  (85) 

FROST  PENETRATION 

We  will  seek  solutions  with  a  positive  V0  in  Sp.  If  such  solutions  exist,  by  eq  P7  P'  (8  -) 
vanishes  and  eq  P2  is  reduced  to 


/(8-)  =  -K2(7i)T'(  5-).  (86) 

Using  eq  PI,  43  and  86,  and  neglecting  sensible  heat  terms,  we  obtain 

e0Pio  (l  -  ViWo1)^  =y~fo  (87) 


where 


i 

e0  —  s2  (1 — r|_1S3 y),  s3  =  koL  /  (didj1  - 1)  (88) 

Vi=v(Ti),  Y -ya0  / r\,  y  =  K2{T1).  (89) 


We  will  write  eq  P4  as 


pa0(l-V1Wo1)Vo=Q-/o  (9°) 

where 

Q  =  (fc1a1-“fc0a0)/L.  (91) 

It  should  be  noted  that  e0  and  Y  are  functions  of  T\  because  y  is  a  function  of 
Using  eq  87  and  90,  we  will  express  V0  and  fQ  in  terms  of  Y  and  Q  as 

eiPlo(l-vW)V0=Q-Y  (92) 

eifo=r-e0Q  (93) 

where  e2  is  a  positive  function  of  Tj  defined  as 
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(94) 


e1  =  l-e0=di1d2[l+(k0f\)-1Ly]. 

Now  the  problem  of  finding  a  solution  with  positive  is  reduced  to  that  of  finding  7\  <  0 

that  satisfies  eq  47, 92  and  93. 

It  follows  from  eq  92  and  93  that  there  are  two  possible  types  of  solutions  satisfying  one 
of  the  following  conditions  given  as 

e0  <  0  or  e0  >  0  and  Y  >  e0Q,  then  f0>  0  (95) 

e0  >  0  and  Y  <  e0Q,  then  fa  <  0.  (96) 

Since  s3  is  a  positive  number  and  K2  (T)  is  an  increasing  and  continuous  function  if  a  =  0, 
we  will  define  Tx  <  0  as 


^2  (Tx )  -  T)  /  s3 ,  a  =  0. 

(97) 

Using  eq  92  and  93,  we  will  write  eq  47  as 

g1(T1)  =  g2(T1) 

(98) 

where  gi  and  g2  are  defined  as 

gi(Tl)  =  a  +  (do/K0)f0(T1) 

(99) 

(100) 

and  Y  and/ are  given  as 

T'(srT i)  =  -ri-i[a0  +  s2s3(i(s,T1){Q  -  Y(Ta)}  /  e^T,)] 

(101) 

f{s,Ti)  =/o(Ta)  +  s2p(s,T1){Q  -  Y(Ta)}  /  e1(T1) 

(102) 

/0(Ti)  =  [Y(Ta)-eo(Ta)Q]/ei(T1) 

(103) 

p(s,Ta)  =  [if0  -  v(s) ]/[w0  -  v(Ti)]. 

(104) 

It  is  noted  that  gi  and  g2  may  be  discontinuous  at  Ti  =  T0  due  to  the  singularity  of  K2  (Ti). 
Below  we  will  study  the  properties  of  gi  (7\)  and  g2  (T{)  that  will  be  used  later  for  existence 
proofs. 


Proposition  4 

For  given  a0  and  oq,  gi  (Tj)  and  g2  (Tj)  have  the  following  properties: 


(Tm)<g2(Tm) 

(105) 

gi(T\)>0  for  Tm  <Ta  <0 

(106) 

g2 (Ti)  <  0  for  Tm  <  Ti  <  0  if  e0  >  0  and  V^,  >  0 

(107) 

where  a  dot  denotes  differentiation  with  respect  to  Tj. 
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Proof 

When  T\  =  Tm,  then  Y  =  0,  e0  =  s2  and  e\  =  1  -  S2-  From  eq  101  through  103,  we  obtain 


fa  ~  ~  s2  Q/  (1  “  sf)  <  0 

(108) 

f/T'  =  s2T|(1  -  p)Q  /  (eaa0  +  s2s3pQ)  >  0. 

(109) 

It  follows  from  eq  108  that  gi  (Tm)  <  0. 

Since  Tm  <  Ta  by  eq  40,  we  will  write  eq  100  as 

it  \  j.  flcr  K-2{s)  ,  f0  f(s,  Tj)  j 

22 1 [Tm  )■ -  a  JTm  Kl(s)  ds  +  lTm  Kl(s)r(S/  Tl)  ds- 

(110) 

Since  the  second  and  the  third  terms  in  the  right  side  of  eq  110  are  positive, 
(Tm)  >  a;  eq  105  holds  true. 

Differentiating  eq  99  with  respect  to  T2,  we  obtain 

we  find  that  g2 

gi(Ti)  =  (?>o/K0)(aQ/r\){ei/ei)y 

(111) 

where  e3  is  defined  as 

e3=e1  +  s2s3(Q-Y)/a0  =  [(l-s2)a0  +  s2s3Q]/a0  >0. 

(112) 

It  follows  from  eq  111  that  106  holds  true. 

Differentiating  eq  100  with  respect  to  Tj,  we  obtain 

22 (Ti)  (il / a0)(e3  / ei  Ka(s)KlH-e4jx(s,T1)J2  ds 

(113) 

where  e4  and  W  are  defined  as 

e4(Ti)  =  (s2s3  /  a0)(Q-Y)/  e1 

(114) 

W(s,Tx)  =  [a0(l-s2|i)+s2s3pQ]Y  +  e0a0ji(Q  -  Y). 

(115) 

Since  V0  >  0,  by  eq  92,  Q  -  Y  >  0.  Differentiating  eq  104  with  respect  to  Tv  we  obtain 


|i(S,T)  =  K  -  vOi)  rV  vOi)  =  Po  (Ti)n  >  o.  (116) 

Since  e0  >  0,  so  W  (s,T\)  >  0.  Hence,  eq  107  holds  true.  □ 

It  is  noted  that  gi(Tj)  and  g2(7i)  maybe  discontinuous  at  T\  =  Ta  due  to  the  singularity 
of  K2  (Tj).  We  will  consider  a  straight  line  L  (<Xis,  a0)  in  Sp  defined  as 

LCais,  «o)  =  [(«•!,  a0):  aa  >  als]  (117) 

where  (ocis,a0)  e  Cs.  We  will  study  the  behavior  of  and  V0  on  L(ais,a0).  Differentiating 

eq  92  with  respect  to  Q,  we  obtain 


ei(io0-Vi)p3oVr0=l-e2fi 


(118) 


where  a  tilde  denotes  differentiation  with  respect  to  Q  for  a  given  a0  and  e2  is  defined  as 
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e2(Ti)  =  (e3M)Y-ji0(Q-Y). 


(119) 


Differentiating  eq  98  with  respect  to  Q,  we  obtain 


£a(Ta)Ti  =£2^) 


(120) 


where  Ex  and  E2  are  defined  as: 


£l(T1)  =  gi(T1)-g2(Ta) 


(121) 


E2(Ti)  =  (e0/e1) 


(80/K0)  +  (ii/a0)[° 


l-R(s,Tj) 


^(s^l  +  ^Ta)]2 


-ds 


(122) 


Now  we  will  begin  our  search  of  solutions  of  eq  92, 93  and  98  in  Sp  with  a  special  case  in 
which  e0  vanishes.  For  the  sake  of  brevity  we  will  refer  the  problem  of  eq  92, 93  and  98  to  as 
Problem  P  hereafter. 


Proposition  5 

There  exists  a  unique  solution  of  Problem  P  such  that  eQ  =  0  and  V0  >  0  if  the  following 
condition  holds  true: 


„  .  f°  K2(s)-(ti/s3) 
hx  K,(s) 


ds. 


(123) 


Proof 

When  e0  =  0,  then  Ti  =  Tx  and  from  eq  93  we  obtain 

/0=Y  =  oc0/s3.  (124) 

Using  eq  101, 102  and  124,  we  will  reduce  eq  98  to 

a0So/(s3Ko)  =  Jj°  S3^s~q-  (125) 

Since  eq  123  holds  true,  for  a  given  a  there  exists  a  unique  and  positive  a0  (a)  that  satisfies 
eq  125. 

We  will  denote  a0  (a)  by  aoc ,  and  consider  a  line  Lc  (ais,aoc)  defined  as 


Mai^oCocHfta^aoc)  :  ai>ais}.  (126) 

Our  aim  is  to  show  V0  >  0  on  Lc.  From  eq  122  we  find  that  E2  (Ta)  vanishes  if  and  only  if  Ti 
=  Tx  or  e0  =  0.  From  eq  111  and  113  we  find  that  gi{Tx)  >  0  and  £2(TX)  <  0;  hence £j  (Tx)  >  0. 
Therefore,  T)  vanishes  in  this  case.  Using  eq  118,  we  find  that  V0  is  positive.  Therefore,  VQ 
is  positive  on  Lc  except  at  a  point  (ais,aoc)  where  yo  vanishes.  O 
We  will  define  g0  (a)  as 


*2(S)-(T1/S3) 

Ki(s) 


ds-  a 


or  we  may  write  g„  as 


(127) 
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s°{°>£mis-MsKmds 


(128) 


where  T0  =  -  o/y  as  defined  by  eq  37.  We  will  define  <rx  as 


<*x  =  -YTx. 


(129) 


It  follows  from  eq  128  that  Tx  <  Ta  or  cx  >  a  when  eq  123  holds  true. 

When  Tx<Ta,  g0  (a)  is  a  decreasing  and  continuous  function  of  a.  We  will  define  ac  as 


?”<a')=fll)‘fe"(,i/S3C  w/s=0‘ or  0t  ■£  K2(t£/S3)*  m 


Ki(s) 


where  Tc  is  defined  as 


Tc  =  -Cc/y>Tx.  (131) 

It  is  noted  that  Tc  is  uniquely  determined  by  the  properties  of  a  given  soil  and  that  Tc  >  Tx  or 
<7C  <  ax.  Using  gor  when  a  <  oc,  we  will  write  aoc  of  Proposition  5  as 


a0C=(s3K0/8o)(cc-c).  (132) 

When  ct  <  gc,  the  Lc  (als,  aoc)  divides  Sp  into  two  regions  Sp  (e0  >  0)  and  Sp  (e0  <  0)  (Fig. 
2b).  The  region  Sp  disappears  if  a  >  crc.  In  terms  of  a0  we  may  define  Sp  and  Sp  as 

SP  ={(ai,a0)6Sp  :  a0>a0C,  eo>0}  (133) 

Sp  {(a!,a0)€Sp  .  (Xq  <cx0C/  c0  <0^.  (134) 

SOLUTIONS  IN  Sp  WHEN  a  <  ac 
We  will  consider  a  line  Lq  defined  as 

Lc(ais,a0)  =  {(a1,a0):aa>als  and  ao>aoc>0}.  (135) 

It  is  clear  that  L+c  belongs  to  Sp  where  Ti  <  Tx. 

Proposition  6 

In  Sp  Ti  <  Tx  <  Ta  if  a  <  ac . 

Proof 

Suppose  that  there  exists  a  point  in  Sp  such  that  Tx  >  Tc.  We  will  write  eq  128  as 

ac=c- y(Tx -T„)-(ri/ s3)^°  ds.  (136) 

It  follows  from  eq  136  that  G  >  gc.  This  contradicts  the  assumption.  □ 

Proposition  7 

Suppose  that  a  solution  of  Problem  P  exists  on  L+c,  then  V0  >  0  and  %  >  0  on  L^. 


17 


Proof 

First  we  will  examine  the  behavior  of  V0  in  a  neighborhood  of  =  ais  on  L+c.  When  ai 
approaches  als/  then  (Q  -  Y)  approaches  zero.  Since  e3  approaches  e\,  eq  118  is  reduced  to: 


ei(ro0-Vi)pioV0  =  l-YT1. 

(13  7) 

Also  eq  120  is  reduced  to 

YT1E3(T1)  =  E4(T1) 

(138) 

where  £3  and  £4  are  defined  as 

(139) 

E,(Ti)  =  e0  l  +  (fC0/50)(ri/a0)^  . 

(140) 

Since  0  <  eQ  <  1,  so  £3  >  £4  >  0  or  YT3  <e0  <1.  Hence  we  find  that  VQ  is  positive  in  a 
neighborhood  of  04  =  ocis.  Suppose  that  there  is  a  point  on  L*  where  V0  <  0.  Since  V0  (Q)  is 
continuous  on  Lc  for  T3  <  Tc,  one  can  find  a  point  on  such  that  V0  vanishes.  However, 
this  contradicts  Prop.  1.  Therefore,  V0  must  be  positive  on  L*.  Since  eo>0  and  V0>Oon  L+, 
in  eq  120  £2  (Tj)  >  0  and  £2  (Ti)  >  0  by  eq  101.  Therefore,  %  >  0  on  (ais ,  a0)  and  on  Cs.  □ 

When  e0  >  0,  a  solution  with  negative/,,  may  exist  by  eq  96.  We  will  study  such  a  possi¬ 
bility  below. 

Proposition  8 

In  Sp  f0  is  positive  if  a  <  ac. 

Proof 

Suppose  that  there  exists  a  solution  Tj  =  Tp  of  Problem  P  on  L+c  (ais ,  aQ)  such  that/0 
(Tp)  =  0  or  Y  (Tp)  ss  eQ  (Tp)Q.  From  eq  101  and  102  we  obtain 


T'(s,Tp )  =  -(1/  T!)[a0  +  s2s3p(s,Tp)Q] 

(141) 

/(s,Tp)  =  s2p(s,Tp)Q. 

(142) 

Tp  is  a  solution  of  the  following  equation  given  as 

a  =  gi{Ti:/0  =  0}  =  g3(Tj) 

(143) 

where  g3  is  defined  as 

on  (Ti )  -  f °  ^  ds  (ti  /  s3 )  f °  ^4P(S/  Ti )  ^ 

hih(s)  ^ ^  3  "Ti  K\ (s)[l  +  e4p(s,  Tiq  S 

(144) 

where  £4  is  defined  as 

e4  =  (s2s3  /  a0)  Q  =  s2s3Y  /  (a0e0). 

(145) 

It  is  easy  to  find 
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ft(T")=<’+£t$'fe>a 


(146) 


Differentiating  eq  144  with  respect  to  7\,  we  obtain 

It  follows  from  eq  147  that  g3 (Ti )  <  0 .  Since  gj(Ti)  is  continuous  and  Tm  < 
a  solution  Tj  =  Tp  of  eq  143  to  exist  we  must  have 


(147) 


Ti  <  Tx  <  Tc,  for 


<r-g3(Tx)>0.  (148) 

When  Tj  approaches  Tx,  e0  approaches  zero  and  g3  (T2)  approaches  <rc.  Hence  eq  148  is 
reduced  to 

a  >  <7C.  (149) 

It  is  clear  that  eq  149  does  not  hold  in  this  case.  Hence  the  solution  Tp  does  not  exist  on  L^. 
Since/0  >  0  in  a  neighborhood  of  oci  =  (Xis  on  L*  and/0  (Q)  is  continuous  on  for  Tj  <  Ta, 
f0  is  positive  on  Lc .  Since  a0  is  an  arbitrary  positive  number,/0  is  positive  in  Sp  if  <j  <  ac.  □ 

Proposition  9 

There  exists  a  unique  solution  of  Problem  P  in  Sp . 

Proof 

We  will  consider  a  line  Lc  defined  by  eq  135.  In  view  of  Proposition  4  we  need  to  show 
that  g2  (Tx)  >  g2  (Tx).  From  eq  99  and  100  we  obtain 


8\(Tx)  =  a05o/(s3K0) 

(150) 

82  (lx)  =  ac- 

(151) 

For  the  sake  of  convenience  we  will  define  WQ  (T )  for  T  <  0  as 

W0(T)=a(T)-a(T). 

(152) 

From  eq  150  and  151  we  obtain 

W0(TX)  =  oc050  /  (s3K0)  +  o,“0’c. 

(153) 

Using  eq  132,  we  will  reduce  eq  153  to 

Wo(Tx)  =  8o(ao-aoc)/(s3Ko)* 

(154) 

Since  a0  >  aoc,  we  find  that  gi  (Tx)  >  g2  (Tx).  Since  gi  (T{)  and  g2  (Ti)  are  continuous  for  Ti  < 
Ta,  there  exists  a  unique  solution  of  Problem  P  on  L+.  Since  a0  is  arbitrary,  for  any  given 
point  (a!  a0)  in  Sp  there  exists  a  unique  solution.  □ 

When  V0  >  0,  unfrozen  water  may  exist  in  R2  (T  <  T{)  and  the  amount  of  unfrozen  water 
in  R2  depends  on  T.  Hence,  the  rate  of  heave  depends  on  T  and  is  given  for  T  <  Ti  as 
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(155) 


d2r{T)  -f0  +  s2  [pio  -  p30v(T)]^o- 
Using  eq  87,  we  will  reduce  eq  155  to 

d2r(T)  =  Y  +  s2p3o[v!  -v(T)  +  (s3/ti)i/v1]  V0  (156) 

where  v3  is  defined  as 

vifli  )  =  w0-v(T1).  (157) 

It  follows  from  eq  156  that  r  is  positive  for  T  >  Tm.  In  engineering  practices  the  frost  heave 
ratio  h  and  the  water  intake  ratio  hw  are  often  used.  These  are  defined  as 

h  =  r(T)/V0  (158) 

K=fo/V0.  (159) 

According  to  M\,  h  and  hw  are  given  as 

h  =  a0y/  (d2r\V0)  +  (s2p30  /  d2)[va  -  v(T)  +  (s3  /  ri)yvx]  (160) 

K=u0y /(r\V0)-e0p3Qvl.  (161) 


SOLUTIONS  IN  Sp  WHEN  <T>oc 

When  cr  >  oc,  the  region  Sp  defined  by  eq  130  disappears,  so  Sp  =  Sp .  We  will  consider  a 
line  L+  (Fig.  2c)  defined  as 

L+(ais,  a0)  =  [(eq,  a0):  aa  >  als  and  a0  >  0]  (162) 

where  (ais,  a0)  e  Cs. 


Proposition  10 

Suppose  that  there  exists  a  solution  T;i  of  Problem  P  on  L+  such  that  V0  >  0.  Then  T-t  <  Ta 
if/o  >  0,  while  T1  may  be  greater  than  Ta  if  f0  <  0. 


Proof 

We  will  write  eq  98  as 


(80  /K0)- 


f° _ ! _ 

Jii  Ki(s)T'(s,  Tj) 


=fti‘,s+fe/e')(C-Y)Jr" 


_ML  * 

K^T'isJi)  ' 


(163) 


Since  V0  >  0,  so  (Q  -  Y)  >  0.  It  follows  from  eq  163  that  T1  <  Ta  if  fQ  >  0,  while  Tj  may  be 
greater  than  Ta  if  f0  <  0.  O 

Since  e0  >0,T1<  Tx.  There  are  two  cases:  Case  1.  ac<  o  <  ox  or  Ta>Tx>  and  Case  2. 
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cr  >  crx  or  Ta  <  Tx.  First  we  will  study  Case  1.  In  Case  1  if  a  solution  of  Problem  P  exists  on 
L+,  then  V0  >  0  and  Tj  >  0  on  L+  by  the  same  reasoning  as  used  in  the  proof  of  Proposition  7. 

Proposition  11 

There  exists  a  unique  solution  <  Tx  of  Problem  P  on  L+  if  c  <  <jx. 

Proof 

For  a  given  point  (oq,  a0)  on  L+,  we  will  search  a  solution  of  Problem  P  under  the  as¬ 
sumption  that  yo>0at  this  point.  Because  of  eq  105, 106  and  101  for  a  unique  solution  Tj  to 
exist  we  must  have 

W0  (Tx)  >  0.  (164) 

Using  eq  99  and  100  we  obtain 


W0(TX)  =  a0S0  /  (s3K0)  +  g-oc  >  0.  (165) 

Hence,  there  exists  a  unique  solution  <  Tx  of  Problem  P  on  L+  if  V0  >  0.  But  if  such  a 
solution  exists,  then  V0  >  0.  Therefore,  there  exists  a  unique  solution  of  Problem  P.  □ 

Proposition  12 

There  exists  a  unique  solution  Ti  =  Tp  of  Problem  P  on  L+  such  that/0  (Tp)  =  0  if  a  <  ax. 

Proof 

It  is  easy  to  see  that  Tp  is  also  a  solution  of  eq  143  of  Prop.  8.  In  this  case  a  >  oc.  Hence, 
there  exists  a  unique  solution  Tp  of  eq  143.  Since/0  =  0,  V0  >  0  by  eq  90.  We  will  show  that 
this  solution  is  located  on  L+.  We  will  denote  the  point  such  that  T\  =  Tp  by  (aip,  a0)  and  the 
point  on  Cs  by  (ais,  a0).  Since  f0  =  0  or  Y  -  e0Q  =  0  at  (aip,  a0),  we  obtain 


aip=(l/fci)[fc0+^2(Tp)/{Tic0(Tp)}]a0  .  (166) 

where  e0  (Tp)  is  given  as 

eo(Tp)  =  s2[l-(s3/r|)K2(Tp)].  (167) 

From  eq  69  we  obtain 

ais  =(l/h)[k0+LK2(Ts)/r\]a0.  (168) 

Using  eq  166  and  168,  we  obtain 

«iP  - «is  =  [K2(Tp) - e0(Tp)i<C2(Ts)]  L  /  [k^e^)].  (169) 


Since  Tx  >0  on  L+,  K2(Tp)  >  K2( Ts).  Also  1>  eo(Tp)>0.  Hence,  alp  >  als.  Therefore,  the 
point  (aip,  a0)  is  on  L+.  □ 

The  unique  solution  Ti  of  Proposition  12  satisfies  eq  143.  We  may  write  eq  143  and  145  as 

(170) 


CT  =  g3  (Tp) 

h(Tp)  =  s2s3K2(T?)/{r\e0(Tv)}. 
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(171) 


It  is  easy  to  see  from  eq  144  that  Tp  depends  on  neither  a0  nor  oq.  Differentiating  eq  170 
with  respect  to  a,  we  find  that  Tp  is  a  decreasing  function  of  a. 

We  will  introduce  a  straight  line  Lp  starting  from  the  origin  (Fig.  2c)  defined  as 


Lp  =|(a1,a0):a0  =fci[fc0+LK2(Tp)/{iie0(Tp)}]"1a1 


(172) 


The  line  Lp  divides  the  region  Sp  into  two  regions  Sx  and  Spp  defined  as 


Sx={(a1,a0):a0<fc1[fc0+LX2(Tp)/{Tic0(Tp)}]-1a1}  (173) 

Spp  ={(ai, a0):a0  >h[k0  +LK2(Tp)/{r1e0(Tp)}]-1a1}.  (174) 


It  is  clear  that/0  <  0  in  Sx  while/0  >  0  in  Spp 

Now  we  will  study  Case  2  where  Ta  <  Tx.  In  this  case  a  solution  7\  may  be  greater  than 
Ta  by  Proposition  10.  First  we  will  search  a  solution  Ti  of  Problem  P  such  that  Ta  <  T0  <  Tx 
on  L+  under  the  assumption  of  V0  >  0.  If  such  a  solution  exists,  then  V0  >  0  and  %  >  0  on  L+ 
by  the  same  reasoning  as  used  in  the  proof  of  Proposition  7. 


Proposition  13 

There  exists  a  unique  solution  of  Problem  P  on  L+  such  that  Tj  <  Tc  <  Tx  and/0  >  0  if 
cr>ox. 


Proof 

First  we  assume  that  V0  >  0.  Because  of  eq  105, 106,  and  107  for  a  unique  solution  Ti  to 
exist  we  must  have: 


Wo(To-)  =  gi(Tc-)-g2(Tc-)>  0.  (175) 

From  eq  99  and  100  we  obtain: 

Wo  - MT ’.-Vo  -fa  /  elm -  Y)  ll  &  (176) 

where  a0  is  defined  as 

»o(l')=(8o/Ko)7(n/a0)gKi(s)[1+1g4)i(tT)1&  (177) 

where  f0,  e1  and  Y  are  functions  of  TG~,  and  e4  defined  by  eq  114  is  given  as 


e4(Tcr“) -(S2S3  /  a0)[Q- Y(Ta-)]/ ei(Tc-).  (178) 

From  eq  177  we  find  that  >  0.  Since  VQ  >  0,  so  (Q  -  Y)  >  0.  When/0  >  0,  Wc  >  0.  Therefore, 
a  unique  solution  T!  of  Problem  P  such  that/0  >  0  exists  on  L+  if  V0  >  0. 

We  will  denote  one  of  such  unique  solution  on  L+  such  that  f0  =  0  by  Tp.  Let  Tp  be 
located  at  (alp,  a0)  on  L+.  Since/0  >  0  in  a  neighborhood  of  Cs  in  Sp ,  the  point  (alp,  a0)  must 
be  in  Sp .  We  will  consider  a  segment  £p  (a0)  of  L+  defined  as 

^P(a0)  =  [(«i/«o):«is  <«i  <«iP]-  (179) 

Since /0  >  0  on  fp(a0),  there  exists  a  unique  solution  T1  <  TCT  on  fp(a0).  Since  VQ  >  0  in  a 
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neighborhood  of  Cs  in  Sp  and  V0  (Q)  is  continuous  on  £p  (a0 ) ,  V0  >  0  on  £p  (a0) .  Therefore, 
there  exists  a  unique  solution  of  Problem  P  on  lp  (a0 )  such  that /0  >  0.  When/0  =  0,  V0  >  0 
because  from  eq  92  and  93  we  find 


Q-Y  =  ei(Q-/o)  =  eiQ>0.  □ 


(180) 


We  will  define  Lp  by  eq  172  where  Tp  is  the  unique  solution  in  Prop.  13.  Then  the  line  Lp 
divides  the  region  Sp  into  two  regions,  Sx  and  Spp  defined  by  eq  173  and  174,  respectively. 
Proposition  13  implies  that  there  exists  a  unique  solution  in  Spp  such  that/0  >  0.  We  will 
search  a  solution  Ti  of  Problem  P  in  Sx  such  that  T2  <  Ta  <  Tx  and/0  <  0.  It  is  easy  to  find 
from  eq  180  that  V0  >  0  when/0  <  0. 

When/0  <  0,  we  will  study  the  behavior  of  W0  (T)  for  Tp  <  T  <  Tx.  When  T  =  Tp,/0  =  0. 
From  eq  176  we  obtain 


(181) 


Using  eq  121,  we  obtain: 


W0(T)  =  gi(T)-g2(T)  =  E1(T).  (182) 

Because  of  eq  106  and  107  we  find 

W0(T)  >  0  for  Tp  <  T  <  Tx .  (183) 

We  may  write  eq  183  as 


W0(Tct-)>  0  for  Tp  <  Ta  <  Tx .  (184) 

We  will  consider  L+  defined  by  eq  162.  We  have  found  that  there  exists  a  unique  solution 
T1  of  Problem  P  for  als  <  ax  <  alp.  As  oq  increases  from  als  to  oqp,  increases  from  Ts  to  Tp 
and/0  decreases  from/s  to  zero.  Because  of  eq  184  there  exists  a  unique  solution  T1  of  Prob. 
P  for  alp  <  cci  such  that  Tp  <  Ti  <  Ta  <  Tx  and/0  <  0.  We  will  state  our  findings  below. 

Proposition  14 

There  exists  a  unique  solution  7^  of  Problem  P  on  L+  with  alp  <  oq  such  that  Tp  <  <  Ta 

and  f0  <  0  if  a  >  ax. 

As  stated  above  there  may  exist  a  solution  Ti  of  Problem  P  such  that  >  T0  in  Case  2. 
But  we  are  not  certain  of  the  existence  of  such  a  solution. 


SOLUTION  IN  S" 

We  will  consider  a  line  Uc  defined  as 

Lc(ais,a0)  =  {(ai,a0)  :ai  >ais  and  aoc>ao>0}  (185) 

where  («iSi  a0)  e  Cs  and  (ai,  aoc)  is  on  the  line  Lc  (Fig.  2b).  It  is  clear  that  L~c  belongs  to  Sp 
where  e0  <  0. 


23 


Proposition  15 

There  exists  at  least  one  solution  Ta  of  Problem  P  on  L~c  such  that  Tx  <  T2  <  Ts  <  Ta,  V0  >  0 
and/0  >  0. 

Proof 

Since  e0  <  0  in  this  case,  we  must  search  a  solution  T:  on  Uc  such  that  T-y  >  Tx.  First  we 
examine  eq  98  when  T2  =  Tx.  Since  a0  <  aoc  in  this  case,  using  eq  154,  we  find  that  g2  (Tx)  < 

&(TX). 

We  will  assume  that  V0  >  0  and/D  >  0  for  the  time  being.  From  eq  99  and  100  we  obtain: 


W0 (Ts)  =  a0 (Ts)/0  -  ai(Ts)  +  (rj/a0 )a2 (Ts) 
where  a0  is  defined  by  eq  177,  and  a2  and  a2  are  defined  as 

fo  |t(s,T) 


Jr  Kl(s) 

a2m-fe/t,)(Q-r)j;KiWp;t;fer)fc 


Since  Ts  is  the  solution  of  eq  98  when  V0  =  0  (Q  =  Y),  we  obtain 

rO  ds 


W0  (%)  = 


(50/K0+(Ti/a0)L 


/s-«i(Ts)  =  0  ifFo=0 


Jrs  Kx(s)JJ 

where /s  is  given  by  eq  68.  Using  eq  189,  we  will  reduce  eq  186  to 


(186) 

(187) 

(188) 

(189) 


W0(TS)  =  (-1  /  e1)(80  /  K0)£4(Ts)(Q  -  Y)  (190) 

where  E4  is  defined  by  eq  140.  It  follows  from  eq  190  that  gy  (Ts)  >  g2  (Ts).  Since  and  g2  are 
continuous  functions  of  Tj  for  T2  <  Ta,  there  exists  at  least  one  solution  T2  of  Prob.  P  on  Uc 
such  that  Tx  <  Tj  <  Ts  <  TCT  if  170  >  0  and/D  >  0. 

When  such  a  solution  exists,  from  eq  137  and  138  we  find  that  %<Q  and  V0  >  0  in  a 
neighborhood  of  oq  =  ocis  in  Sp .  By  the  similar  reasoning  to  that  used  in  the  proof  of  Prop¬ 
osition  7,  we  find  that  VQ  >  0  on  L7 .  Since  V0  >  0  and  Ty  <  Tc  on  L~ ,  hence  by  Proposition  10 

fa  >  0.  □ 

It  is  clear  that  a  solution  of  Proposition  15  is  unique  if  gj  (T\ )  <  0  on  Lc .  We  will  study 
the  sign  of  g2  (Ty )  below.  We  will  write  eq  113  as 


^(Ti)  =  -(l/a0)(e3/e12)W1(T1) 
where  Wj  is  defined  as 


Wi(Ta)  = 


.(n/«o,li^ 


■ds. 


i(s)(l  +  e4p)2 

The  sign  of  g2  (Tj )  depends  on  that  of  Wi(Ti).  We  will  write  eq  192  as 


(191) 


(192) 


Wi(Ti)  =  Wi<x0  +  w2  a4. 


(193) 


w\  and  w2  are  defined  as 


(fl)  =  W3  +  (ele5  -  y)a 4  >  0 


(194) 
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w2  (Ti  )  =  Si(y-e5  )a4  =  [s^  /  (s2s3v1)]e6a4 

(195) 

si  =(h/ko)(l-s2)>0 

(196) 

^5  (?) ) —  — ^lsop0  /  (s2s3)  >  0 

(197) 

e6(T1)  =  ^-(-e0v1) 

(198) 

fl3(Tl)-JTlKl(s)(1  +  e4[t)2ds 

(199) 

(200) 

The  sign  of  w2  is  the  same  as  e6  and  eg  is  a  property  of  a 
the  case  where  the  following  condition  holds  true: 

given  soil.  First  we  will  consider 

e6(Ta)>0  for  T4<TS. 

(201) 

Suppose  that  T,  is  a  solution  of  Proposition  15.  Then  W1(T1)  >  0  and  g2 (Tj)  <  0  on  L~  in  this 
case.  Therefore,  this  solution  is  unique.  From  eq  120  we  find  that  T4  <  0  on  Vc . 

Next  we  will  study  the  case  in  which  eq  201  does  not  hold  true.  We  will  examine  the 
behavior  of  W1(T1)  on  Lc .  When  cxj  approaches  alg,  approaches  Ts,  als  is  given  by  eq  168. 
Using  eq  88,  we  reduce  eq  168  to 


«is  —  (eis/si)a0  (202) 

where 

ei»  =  ei  (Ts).  (203) 


We  will  write  Wj  as 

Wi(Ti)  =  («3  -  e0a4)ya0  +  (s1a1  -  eia0)(y  -  e5)a4  (204) 

when  (Xj  approaches  als,  the  second  term  of  the  right-hand  side  of  eq  204  vanishes  and 
Wj(Tj)  approaches  W^T)  given  as 

Wi(Ts)  =  (a3  -  e0a4)ya0  >  0.  (205) 

It  follows  from  eq  205  that  W^Tj)  >  0  in  a  neighborhood  of  =  als  on  Lc.  It  is  easy  to  find 
that  Wx  vanishes  when  a:  becomes  infinite. 

Since  the  second  term  of  the  right-hand  side  of  eq  204  is  negative  on  L~,  it  is  possible 
that  W1(T1)  may  become  negative  at  some  point  on  L~.  If  such  is  the  case,  then  there  ex¬ 
ists  at  least  one  point  (alg,  ac)  on  Lc  where  Wa  vanishes  because  is  a  continuous  func¬ 
tion  of  ar  From  eq  193  we  obtain 


aig  =-(a?i  /w2)a0. 


(206) 


From  eq  202  and  206  we  obtain 


(dig  -  ais)  /  a0  =  a3Tie7[l-  (a4  /  a3)eos]  /  (— s2s3vjw2) 


(20 7) 
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where  e7  is  defined  as 


e7(Ti)  =  e6  -e0v1[l-(a4  / a3)e0\  /[l-(a4  /  a3)eos]  (208) 

e0s=e0(Ts).  (209) 

It  follows  from  eq  207  that  alg  >  als  if  e7  is  positive.  This  implies  that  there  exists  a  point 
(dig,  a0)  on  L~c  such  that  Wj  vanishes  and  that  Wi  maybe  negative  for  04  >  alg.  When  Wj  is 
negative,  then  g20\)  is  positive  and  the  uniqueness  of  the  solution  is  not  warranted. 
Hence,  if  e7  is  positive,  the  solution  of  Proposition  15  is  unique  under  condition  given  as 

eg(Ti)  <  0  and  e7(Ti)  >  0  for  Ti  <  Ts  and  ai  <  alg.  (210) 

On  the  other  hand,  if  e7  <  0,  then  alg  <  cci8.  This  clearly  implies  that  Wj  >  0  on  Uc  and  that  a 
solution  of  Proposition  15  is  unique.  We  will  present  our  finding  by  the  following  proposi¬ 
tion. 

Proposition  16 

The  solution  of  Prop.  15  is  unique  if  either  e6  (T{)  >  0  or  e7  (T2)  <  0  for  Tx  <  Ts.  When  e7  (Tj) 
>  0  for  Ti<Ts,  the  solution  is  unique  if  eq  210  holds  true. 


APPLICATIONS 

We  will  describe  the  use  of  traveling  wave  solutions  obtained  above  for  the  empirical 
verification  of  the  model  Mi  below.  It  is  known  (Andersland  and  Anderson  1978)  that  the 
empirically  determined  function  v(T)  under  equilibrium  conditions  takes  a  form  given  as: 


v(T)  =  A0|T|-^i  for  T  <  0 


(211) 


where  A0  and  A\  are  positive  constants.  Experimental  methods  were  proposed  to  deter¬ 
mine  Ki  (Williams  and  Burt  1974,  Horiguchi  and  Miller  1983)  and  K2  (Perfect  and  Williams 
1980).  Horiguchi  and  Miller  (1983)  empirically  found  that  fC]  of  several  frozen  porous  me¬ 
dia  also  takes  the  same  form  as  eq  211.  Since  v  and  Kj  are  known  to  be  bounded,  we  will 
use  forms  given  as 


[w0(A/Tp 


A  <T  <0 

ifVo>0 

A  >  T 


(212) 


Ki(T)  = 


K0 

[K0(A/T)fci 


A  <T  <0 
A  >T 


(213) 


where  A  is  a  small  negative  number,  b\  and  b3  are  positive  numbers.  When  VQ  =  0,  v  (T)  is 
not  needed  in  the  balance  equations  of  mass  and  heat.  However,  K,  (2  =  1,2)  implicitly  takes 
account  of  the  composition  of  the  frozen  fringe. 

Recently,  Nakano  and  Takeda  (1994)  empirically  found  that  K2  (T)  of  Kanto  loam  can  be 
described  in  the  same  form  as  eq  211  for  T  <  Ta.  Using  eq  37,  we  will  describe  K2  as 


{K20=yKo 


K2(T)  =  \ 


y  Ki(T) 
K20(A  /  T)b2 


A  <T  <0 
T0  <  T  <  A 
T  <T0. 


(214) 
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The  empirically  determined  values  of  parameters  in  eq  212, 213,  and  214  for  Kanto  loam 
(Nakano  and  Takeda  1994)  are  wQ  =  0.740,  A  =  -1.5  x  10  ~4  °C,  b3  =  0.110,  Ka  =  1.77  x  103  g/ 
(cm-d-MPa),  bx  =  0.520,  K2o  =  1.98  x  103  g/(cnvd-°C)  and  b2  =  1.04. 

The  existence  of  the  boundary  Cs  (Fig.  2a)  has  been  verified  empirically  for  three  types 
of  soils  including  Kanto  loam  under  o  =  0  (Takeda  and  Nakano  1990).  According  to  Prop.  2, 
Ts  is  a  decreasing  function  of  o.  This  implies  that  the  region  Sj  decreases  as  <j  increases, 
which  is  also  verified  by  the  data  of  Kanto  loam  (Takeda  and  Nakano  1993).  We  will  con¬ 
sider  a  freezing  test  in  which  a  soil  sample  with  a  uniform  initial  temperature  Ta  >  0°C  is 
frozen  from  the  bottom  up  while  the  bottom  temperature  is  kept  constant  at  Tb  <  0°C  and 
the  top  temperature  at  Ta.  The  temperature  field  changes  rapidly  at  the  start  of  the  test. 
However,  as  time  elapses,  the  rate  of  the  change  slows  down  so  that  the  transient  freezing 
may  be  approximated  by  a  series  of  quasi-steady  freezing  steps.  Hence,  the  later  part  of  the 
experiment  can  be  represented  by  a  trajectory  in  Figure  3,  consisting  of  points 
a(f)  =  {ai(f),a0(f)}  for  t2>t>t0. 
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A  point  a  (f0)  is  in  Sp  where  frozen  soil  without  any  visible  ice  layer  grows.  As  oq 
decreases  and  a0  increases  with  time,  the  trajectory  approaches  the  vicinity  of  a(q)  in  Sp. 
The  pattern  of  ice-rich  frozen  soil  grown  in  this  vicinity  evidently  depends  on  the  soil  type 
and  the  magnitude  of  a]  (or  a0)  (Takeda  and  Nakano  1990).  The  results  of  tests  on  Kanto 
loam,  for  instance,  clearly  indicate  that  the  pattern  of  rhythmic  ice  banding  is  formed  at  the 
small  values  of  oq,  while  soil  particles  or  small  aggregates  of  soil  particles  are  evenly  dis¬ 
persed  at  the  greater  values  of  oq.  The  mechanism  of  such  pattern  formation  is  not  well 
understood,  but  possible  causes  include  the  heterogeneity  of  soils  and  instabilities  in  the 
coupled  heat  and  mass  transport  process  with  phase  change.  The  quasi-steady  solution 
described  above  predicts  the  average  ice  content  of  frozen  soil  and  accounts  for  the  growth 
of  an  ice  layer  but  not  the  time-dependent  formation  of  patterns. 

When  a (t)  reaches  a(ti)  on  Cs,  an  ice  layer  emerges.  While  a (t)  moves  toward  a(t2)  on  Le, 
the  growth  of  the  ice  layer  continues  with  decreasing  growth  rate  until  a(t)  reaches  a(t2)  on 
L&  where  the  ice  layer  stops  growing.  At  a(t2)f0  vanishes  and  we  find  from  eq  59  and  60  that 
eq  76  is  reduced  to 


<y  =  -yTa. 


(215) 


27 


This  equation  is  often  called  the  generalized  Clausius-Clapeyron  equation  that  describes  the 
equality  of  chemical  potentials  of  ice  and  water  subjected  to  two  different  pressures.  Radd 
and  Oertle  (1973)  empirically  validated  eq  215. 

The  soil  water  is  known  to  be  expelled  from  the  frost  front  under  certain  conditions 
depending  upon  soil  type,  stress  level,  etc.,  when  a  frost  front  advances  through  a  saturated 
soil.  Such  a  phenomenon  is  often  called  the  pore  water  expulsion,  and  a  concise  review  of 
papers  on  the  subject  was  written  by  McRoberts  and  Morgenstem  (1975).  We  have  found 
that  the  unique  solution  of  Problem  P  exists  such  that  f0  <  0  if  a  >  ac.  It  is  clear  that  the  pore 
water  expulsion  occurs  in  this  solution.  We  will  examine  the  accuracy  of  M\  by  using  exper¬ 
imental  data  of  Kanto  loam  on  water  expulsion  below. 

Takashi  et  al.  (1978)  conducted  numerous  freezing  tests  similar  to  the  test  described 
above  on  overconsolidated  samples  of  silt  and  clay  to  determine  empirical  descriptions  for 
the  heave  ratio  h  and  the  water  intake  ratio  hw.  In  their  tests,  Ta  >  0°C  was  kept  constant  at  a 
value  0.2-0.3°C  higher  than  the  freezing  point  of  the  sample,  while  Tb  (t)  was  decreased  with 
time  from  the  initial  value  Tb  (0)  =  Ta  in  such  a  manner  that  V0  was  kept  nearly  constant. 
After  each  test  h  and  hw  were  determined  by  measured  total  amounts  of  heave  and  water 
intake,  respectively,  for  a  given  set  of  <r  and  V0.  The  empirical  descriptions  obtained  are  giv¬ 
en  as 


h  =  (m1/o)[l  +  (m2/V0)i/2]  +  m0  (216) 

K  =  /o)[l  +  (m2/  Vo)1/2]  -  s2m3  (21 7) 

where  ( i  =  0,1,...,  3)  are  positive  numbers  that  depend  on  a  given  soil.  The  sets  of  con¬ 
stants  nii  for  a  few  kinds  of  soils  have  been  reported  (Ohrai  and  Yamamoto  1991).  Ryokai 
(1985)  determined  the  set  of  constants  m,  for  Kanto  loam  by  a  series  of  freezing  tests  similar 
to  those  of  Takashi  et  al.  (1978).  The  values  of  wt;  are  m0  =  0.0002,  m1  =  0.980  kPa,  m2  =  8.07  x 
103  cm/ d  and  m3  =  0.439.  In  his  tests  (Ryokai  1985)  the  height  of  samples  was  2.0  cm  and  Ta  = 
0.5°C. 

In  the  Takashi's  freezing  test,  a  and  V0  are  constants  but  a0  and  80  vary  with  time.  For 
instance,  in  the  test  by  Ryokai  (1985),  the  value  of  S0  was  2.0  cm  at  the  start  and  decreased 
with  time.  The  value  of  a0  was  0.25°C/ cm  at  the  start,  increased  to  about  1.0°C/ cm  when  a 
quarter  of  a  sample  remained  unfrozen,  and  increased  further  with  time.  It  has  been  recog¬ 
nized  (Ohrai  and  Yamamoto  1991)  that  the  empirical  formulas  (eq  216  and  217)  must  be 
applied  for  cases  where  V0  is  greater  than  about  1.5  cm/d,  because  the  behavior  of  these 
formulas  as  V0  approaches  zero  is  incompatible  with  empirical  findings.  It  follows  from  eq 
216  and  217  that  r  and/0  vanish  as  V0  vanishes.  But  in  reality  when  V0  vanishes,  an  ice  layer 
begins  growing  so  that  r  and/0  do  not  vanish.  It  is  also  important  to  mention  that  the  empir¬ 
ical  formulas  must  be  applied  for  cases  where  a  is  greater  than  about  50  kPa. 

Ideally,  the  results  of  Takashi's  tests  should  be  compared  with  the  theoretical  predictions 
based  on  the  solution  under  the  same  initial  and  boundary  conditions  as  those  of  actual 
tests.  Since  such  unsteady  solutions  are  not  yet  known,  we  will  use  eq  160  and  161  based  on 
the  traveling  wave  solutions  studied  above.  Differentiating  eq  98  with  respect  to  a,  we  ob¬ 
tain 


£i(Ti)^-  =  -l.  (218) 

Since  £i  is  positive,  T3  and  y  are  decreasing  functions  of  cr.  It  follows  from  eq  160  and  216 
that  his  a  positive  and  decreasing  function  of  a  and  V0.  From  eq  161  and  217  we  find  that  hw 
is  also  a  decreasing  function  of  <r  and  V0  and  that  hw  becomes  negative  when  cr  and  V0  be¬ 
come  large,  namely,  water  expulsion  occurs. 
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Assuming  that  r)  =  1,  we  will  calculate  several  important  parameters  of  Kanto  loam  as 

s3  =  4.25  cm  •  d  •  °C/g,  ac  =  453kPa,  <rx  =  998  kPa.  (219) 

We  anticipate  that  water  expulsion  occurs  if  a  >  453  kPa.  In  order  to  calculate  h  and  hy,  by  eq 
160  and  161,  respectively,  we  must  determine  Tx.  Using  eq  212, 213  and  214,  we  will  reduce 
eq  98  to 

F{y(Ti),ao,Vo,c,5o}  =  0.  (220) 

A  detailed  description  of  F  is  given  elsewhere  (Nakano  and  Primicerio  1995).  Since  eq  220  is 
a  nonlinear  algebraic  equation,  for  given  a0 ,  S0,  a  and  V0,  Ta  was  calculated  numerically  by 
the  Newton-Raphson  method. 

Calculating  h  and  as 
functions  of  V0  for  various 
sets  of  (a0/  50,  a)  with  the 
ranges  of  0.1  <  a0  <  1.0°C/ 
cm,  0.5  <  50  <  5.0  cm,  and  0 
<  a  <  1.5  MPa,  we  have 
found  that  the  dependence 
of  h  and  on  a  is  the 
strongest,  and  then  less 
strong  on  Vor  aQ  and  S0  in 
order  of  decreasing  depen¬ 
dence.  The  value  of  80  is 
proportional  to  the  resis¬ 
tance  against  the  flow  of 
water  in  R0.  On  the  other 

hand,  the  flow  resistance  figure  4.  Calculated  values  ofh  (%)  vs.  V0  (cm/d)  under  four  dif- 
of  Ri  increases  with  in-  ferent  values  of  aQ  with  50  =  1.0  cm  and  a  =  500  kPa .  Circles  are 
creasing  a.  When  the  resis-  calculated  by  an  empirical  formula  (Ryokai  1985). 
tance  of  R1  becomes  much 
greater  than  that  of  R0,  the  12 
effect  of  50  diminishes.  The 
effect  of  S0  was  found  neg-  10 
ligible  when  a  is  greater 
than  300  kPa.  The  calculat- 

0 

ed  values  of  h  vs.  V0  and 
//w  vs.  V0  for  Kanto  loam  c? 
under  the  condition  of  a  = 

500  kPa,  and  S0  =  1.0  cm 
with  four  different  values  4 

of  a0  are  presented  in  Fig¬ 
ure  4  and  5,  where  circles  2 
are  values  calculated  by 
the  empirical  formulas  for  0 

Kanto  loam  determined  by 

Ryokai  (1985).  From  Fig-  Figure  5.  Calculated  values  of  hw  (%)  vs.  V0  (cm/d)  under  four 
ure  4  and  5  we  find  that  the  different  values  ofa0  with  50  =  1.0  cm  and  a  =  500  kPa.  Circles  are 
effect  of  ac  on  h  and  hw  is  calculated  by  an  empirical  formula  (Ryokai  1985). 
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Figure  6.  Calculated  values  of  h  (%)  vs.  V0  (cm/d)  with  a0  =  0.75  °C/cm, 
80  =  1.0  cm,  a  =  0.75,  and  1.0  MPa.  Circles  are  calculated  by  an  empiri¬ 
cal  formula  (Ryokai  1985). 


Figure  7.  Calculated  values  of  f0  \g/(cm2  -d) ]  vs.  V0  (cm/d)  with  a0  = 
0.75°C/cm  and  80  =  1.0  cm  and  a-  0.75, 1.0, 1.5 MPa,  and  ».  Circles 
are  calculated  by  an  empirical  formula  (Ryokai  1985). 


significant  and  that  the  calculated  curves  with  a0  =  0.75°C/cm  agree  well  with  the  empiri¬ 
cal  formulas. 

Assuming  that  a0  =  0.75°C/cm  and  80  =  1.0  cm,  we  calculated  y,  h  and  hw  as  functions  of 
Vo  with  a  =  0.75  and  1.0  MPa.  In  Figure  6  predicted  curves  of  h  vs.  V0  with  a  =  0.75  and  1.0 
MPa  are  presented.  From  this  figure  we  find  that  the  predicted  curves  of  h  vs.  V0  tend  to 
deviate  from  the  empirical  formulas  when  Vc  >  6.0  cm/ d.  The  calculated  values  of  f0  vs.  V0 
are  presented  in  Figure  7  when  a0  =  0.75°C/cm,  S0  =  1.0  cm,  and  a  =  0.75, 1.0, 1.5  MPa  and 
°°.  Circles  in  the  figure  are  values  calculated  by  the  empirical  formula.  The  values  of  V0  at/0 
=  0  are  8.94  (7.77)  cm/d,  4.19  (4.31),  and  2.01  (1.88)  for  a  =  0.75, 1.0  and  1.5  MPa,  respec¬ 
tively,  where  numbers  in  parentheses  are  calculated  by  the  empirical  formula.  The  agree¬ 
ment  between  the  predicted  and  empirical  values  of  V0  at/0  =  0  is  satisfactory. 
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